#include "stdafx.h"
#include <cstdio>
#include <cstdlib>  // ham exit()
#include <cmath>
#include <math.h>
#include <string.h>

#include "FFT.h"

FFT::FFT(){};

FFT::FFT(int n)
{
	this->n = n;
	this->m = (int)( log((double)n) / log((double)2) );
	
	// Make sure n is a power of 2
	if(this->n != (1<<this->m))
	{
	}
	
	// precompute tables
	else{
		_cos = new double[n/2];
		_sin = new double[n/2];
			
		for(int i=0; i<n/2; i++) {
			_cos[i] = cos(-2 * PI * i/n);
			_sin[i] = sin(-2 * PI * i/n);
		}
		
		makeWindow();
	}
};

FFT::~FFT() 
{} 


void FFT::makeWindow() 
{
	// Make a blackman window:
	// w(n)=0.42-0.5cos{(2*PI*n)/(N-1)}+0.08cos{(4*PI*n)/(N-1)};
	this->window = new double[this->n];
	//int length = sizeof(window)/sizeof(window[0]);
	for(int i = 0; i < this->n; i++)
		window[i] = 0.42 - 0.5 * cos(2 * PI * i / (this->n-1)) 
					+ 0.08 * cos(4 * PI * i / (this->n-1));
}

void FFT::makeWindow(char* WindowName)
	{
		//Initialize Window
		this->window = new double[this->n];
		int window_type;
		if (strcmp(WindowName, "Hammming"))
			window_type = 1;
		if (strcmp(WindowName, "Hanning"))
			window_type = 2;
		else
			window_type =0;
		switch (window_type)
		{
			//Hamming "Pobre"
		case 1: for (int i = 0; i < this->n; i++)
							this->window[i] = 0.54 - 0.46 * cos(2 * PI * i / (this->n - 1));
			break;

			//Hanning "Good"
		case 2: for (int i = 0; i < this->n; i++)
							this->window[i] = 0.5 - 0.5 * cos(2 * PI * i/(this->n-1));
			break;

			//blackman "Best"
		default: for (int i = 0; i < this->n; i++)
					 this->window[i] = 0.42 - 0.5 *cos(2 * PI * i / (this->n - 1)) +
					 0.08 * cos(4 * PI * i / (this->n - 1));
			break;
		}
	}










double* FFT::getWindow() 
{
	return this->window;
}

void FFT::fft(double* x, double* y)
{
	int i,j,k,n1,n2,a;
	double c,s,t1,t2;					 

	// Bit-reverse
	j = 0;
	n2 = this->n/2;
	for (i=1; i < this->n - 1; i++) {
		n1 = n2;
		while ( j >= n1 ) {
			j = j - n1;
			n1 = n1/2;
		}
		j = j + n1;

		if (i < j) {
			t1 = x[i];
			x[i] = x[j];
			x[j] = t1;
			t1 = y[i];
			y[i] = y[j];
			y[j] = t1;
		}
	}

	// FFT
	n1 = 0;
	n2 = 1;

	for (i=0; i < m; i++) {
		n1 = n2;
		n2 = n2 + n2;
		a = 0;

		for (j=0; j < n1; j++) {
			c = _cos[a];
			s = _sin[a];
			a +=  1 << (m-i-1);

			for (k=j; k < this->n; k=k+n2) {
				t1 = c*x[k+n1] - s*y[k+n1];
				t2 = s*x[k+n1] + c*y[k+n1];
				x[k+n1] = x[k] - t1;
				y[k+n1] = y[k] - t2;
				x[k] = x[k] + t1;
				y[k] = y[k] + t2;
			}
		}
	}
} 

//// Test the FFT to make sure it's working
//#ifdef TEST_FFT
///*
//public static void main(String[] args) {
//	int N = 8;
//
//	FFT fft = new FFT(N);
//
//	double* window = fft.getWindow();
//	double* re = new double[N];
//	double* im = new double[N];
//
//	// Impulse
//	re[0] = 1; im[0] = 0;
//	for(int i=1; i<N; i++)
//	{
//		re[i] = 0;
//		im[i] = 0;
//	}
//	beforeAfter(fft, re, im);
//
//	// Nyquist
//	for(int i=0; i<N; i++) {
//		re[i] = pow(-1, i);
//		im[i] = 0;
//	}
//	beforeAfter(fft, re, im);
//
//	// Single sin
//	for(int i=0; i<N; i++) {
//		re[i] = cos(2*Math.PI*i / N);
//		im[i] = 0;
//	}
//	beforeAfter(fft, re, im);
//
//	// Ramp
//	for(int i=0; i<N; i++) {
//		re[i] = i;
//		im[i] = 0;
//	}
//	beforeAfter(fft, re, im);
//
//	long time = System.currentTimeMillis();
//	double iter = 30000;
//	for(int i=0; i<iter; i++)
//		fft.fft(re,im);
//	time = System.currentTimeMillis() - time;
//	printf("Averaged " + (time/iter) + "ms per iteration\n");
//}
//
//protected static void beforeAfter(FFT fft, double* re, double* im) {
//	printf("Before: \n");
//	printReIm(re, im);
//	fft.fft(re, im);
//	printf("After: \n");
//	printReIm(re, im);
//}
//
//protected static void printReIm(double* re, double* im) {
//	printf("Re: [");
//	for(int i=0; i<re.length; i++)
//		printf(((int)(re[i]*1000)/1000.0) + " ");
//
//	printf("]\nIm: [");
//	for(int i=0; i<im.length; i++)
//		printf(((int)(im[i]*1000)/1000.0) + " ");
//
//	printf("]\n");
//}
//
//*/
//#endif